{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Construct Hamiltonian \n",
    "\n",
    "*Copyright (c) 2021 Institute for Quantum Computing, Baidu Inc. All Rights Reserved.*\n",
    "\n",
    "## Outline\n",
    "\n",
    "This tutorial will demonstrate how to construct a system's Hamiltonian using Quanlse. The outline of this tutorial is as follows:\n",
    "\n",
    "- Introduction\n",
    "- Preparation\n",
    "- Define Hamiltonian\n",
    "- Add control waveforms\n",
    "- Simulation and related tools\n",
    "- Summary"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Introduction\n",
    "\n",
    "Generally, when modeling a superconducting qubits system, we need to define the system's Hamiltonian $\\hat{H}_{\\rm sys}$ firstly. The Hamiltonian, as the total energy, can be used to describe the behavior of entire system:\n",
    "$$\n",
    "\\hat{H}_{\\rm sys} = \\hat{H}_{\\rm drift} + \\hat{H}_{\\rm coup} + \\hat{H}_{\\rm ctrl}.\n",
    "$$\n",
    "It typically contains three terms - the time-independent drift term describing the individual qubits in the system, the coupling term describing the interaction between qubits, and the time-dependent control term describing the control driving acting on qubits.\n",
    "\n",
    "After the Hamiltonian constructed, we can simulate the evolution of the quantum system by solving Schrödinger equation in the Heisenberg picture and then obtain the time-ordered evolution operator $U(t)$,\n",
    "\n",
    "$$\n",
    "i\\hbar\\frac{{\\rm \\partial}U(t)}{{\\rm \\partial}t} = \\hat{H}(t)U(t).\n",
    "$$\n",
    "\n",
    "A variety of functions and a complete set of pre-defined operators are provided in `Quanlse`, which allows users to construct the Hamiltonian for large quantum systems with ease.\n",
    "\n",
    "## Preparation \n",
    "\n",
    "After you installed Quanlse successfully, you can run the program below following this tutorial. Before running this particular program, you would need to import the necessary packages from Quanlse and a constant $\\pi$ from the commonly-used Python library `math`:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "from Quanlse.QHamiltonian import QHamiltonian as QHam\n",
    "from Quanlse.QOperator import number, driveX, driveY, duff\n",
    "from Quanlse.QWaveform import gaussian, dragY1\n",
    "\n",
    "from math import pi"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Define Hamiltonian\n",
    "\n",
    "In this section, we will take a three-level superconducting quantum system consisting of two qubits to demonstrate how to construct a Hamiltonian using Quanlse. Now, we will create such a Hamiltonian in Quanlse by individually adding the time-independent drift terms, the coupling terms and finally the time-dependent control terms.\n",
    "\n",
    "The first step is to instantiate a object 'ham' from class `QHam()`, which contains the basic information regarding the system. Its arguments include the number of qubits `subSysNum`, the system's energy level `sysLevel`, and the sampling-time `dt` (in nanoseconds)."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "ham = QHam(subSysNum=2, sysLevel=3, dt=1.0)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Upon Quanlse 2.0.0, the class `QHam()` allows the users to define different energy levels for qubits in the same system. The users could pass the qubits' energy levels as a list to `sysLevel`. For instance, when `subSysNum=2`, we can pass `[2, 10]` to `sysLevel` to define a system that the qubit 0 has 2 energy levels and qubit 1 has 10. This system is defined in 20-dimensional space $\\mathcal{H}^{2} \\otimes \\mathcal{H}^{10}$ ($\\mathcal{H}^n$ indicates the $n$-dimensional Hilbert space).\n",
    "\n",
    "Then, we add the drift Hamiltonian, which includes the detuning and the anharmonicity terms:\n",
    "$$\n",
    "\\hat{H}_{\\rm drift} = \\sum_{i=0}^1(\\omega_i-\\omega_d)\\hat{a}_i^\\dagger\\hat{a}_i + \\sum_{i=0}^1 \\frac{\\alpha_i}{2} \\hat{a}_i^{\\dagger}\\hat{a}_i^{\\dagger}\\hat{a}_i\\hat{a}_i.\n",
    "$$ \n",
    "\n",
    "These terms are added into `ham` by the function `addDrift()`. We pass in the according parameters and the operator which are pre-defined in the module `QOperator`. ([Clik here](https://quanlse.baidu.com/api/Quanlse/Quanlse.QOperator.html) to view all pre-defined operators)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# The qubit frequency for the qubits, GHz * 2 * pi\n",
    "wq = [4.887 * (2 * pi), 4.562 * (2 * pi)]  \n",
    "\n",
    "# Anharmonicity for the qubits, GHz * 2 * pi\n",
    "anharm = [- 0.317 * (2 * pi), - 0.317 * (2 * pi)]\n",
    "\n",
    "# The drive pulse is in resonance  \n",
    "wd = wq[0]"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Since there are two qubits in this system, we can use a `for` loop of two `addDrift()` functions to add the drift terms. The first parameter is the operator which is a instance of class `QOperator`; `onSubSys` indicates which qubit(s) the term corresponds to; and `coef` is the coefficient for this system; `name` is a user-defined identifier of this term."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Add drift terms to the Hamiltonian\n",
    "for q in range(2):\n",
    "    ham.addDrift(number, onSubSys=q, coef=wq[q] - wd, name='number%d' % q)\n",
    "    ham.addDrift(duff, onSubSys=q, coef=anharm[q] / 2, name='duff%d' % q)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Then, let us add the coupling terms describing the interaction between qubits into system."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Add coupling\n",
    "The coupling term can be written as:\n",
    "$$\n",
    "\\hat{H}_{\\rm coup} = \\frac{g_{01}}{2} (\\hat{a}_0\\hat{a}_1^\\dagger + \\hat{a}_0^\\dagger\\hat{a}_1).\n",
    "$$ \n",
    "\n",
    "In Quanlse, only one line of code is needed to add the coupling term - we use the function `addCoupling()`, which needs to select the qubits' indices `onSubSys` you want and specify a coupling strength `g`. Note that the argument `g` is the coefficient in the front of the coupling term, thus a constant $1/2$ is needed."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "ham.addCoupling(onSubSys=[0, 1], g=0.0277 * (2 * pi) / 2)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Lastly, we add the control terms."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Add control waveforms\n",
    "In a superconducting quantum system, the control Hamiltonian represents external control (microwave, magnetic flux etc.) acting on the qubits. To act X-Y control on qubit 0, the corresponding control terms are given as below:\n",
    "$$\n",
    "\\hat{H}_{\\rm ctrl} = A_0^x(t)\\frac{\\hat{a}_0+\\hat{a}_0^\\dagger}{2} + iA_0^y(t)\\frac{\\hat{a}_0-\\hat{a}_0^\\dagger}{2}. \n",
    "$$ \n",
    "In Quanlse, users can use function `addWave()` to define and add wave to system. Additionally, users can use `clearWaves()` to clear all waveforms in the specified control term.\n",
    "\n",
    "Here, we will take `addWave()` to define and add the control wave. Each waveform function $A(t)$ can be defined using four arguments: start time `t0`, duration `t` and the corresponding parameters `a` 、`tau` 、`sigma`. The function `addWave()` allows us to define waveforms in different ways:\n",
    "\n",
    "- **using preset waveform functions:**\n",
    "Users can use the preset waveforms, here we use the gaussian wave `gaussian` as an example. The supported waveforms are listed here: [API](https://quanlse.baidu.com/api/Quanlse/Quanlse.QWaveform.html).\n",
    "```python\n",
    "ham.addWave(driveX(3), onSubSys=0, waves=gaussian(t0=0, t=20, a=1.1, tau=10, sigma=4), name = 'q0-ctrlx')\n",
    "```\n",
    "\n",
    "- **using user-defined wave functions:**\n",
    "Users can also define a function in the form of `func(_t, args)`, where the first parameter `_t` is the time duration and `args` is the pulse parameters.\n",
    "```python\n",
    "def userWaveform(t0: Union[int, float], t: Union[int, float], a: float, tau: float, sigma: float,\n",
    "             freq: float = None, phi: float = None) -> QWaveform:\n",
    "    \"\"\"\n",
    "    Return a QWaveform object of user-defined wave.\n",
    "    \"\"\"\n",
    "\n",
    "    def func(_t, args):\n",
    "        _a, _tau, _sigma = args\n",
    "        if _sigma == 0:\n",
    "            return 0\n",
    "        pulse = _a * exp(- ((_t - _tau) ** 2 / (2 * _sigma ** 2)))\n",
    "        return pulse\n",
    "\n",
    "    wave = QWaveform(f=func, t0=t0, t=t, args=(a, tau, sigma), freq=freq, phase=phi)\n",
    "    wave.name = \"user-defined wave\"\n",
    "    return wave\n",
    "ham.addWave(driveX(3), onSubSys=0, waves=userWaveform(t0=0, t=20, a=1.1, tau=10, sigma=4), name = 'q0-ctrlx')\n",
    "```\n",
    "\n",
    "In this example, we add the pre-defined Gaussian and DRAG waveforms to the X and Y control terms respectively:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "ham.appendWave(driveX, onSubSys=0, waves=gaussian(t=20, a=1.1, tau=10, sigma=4))\n",
    "ham.appendWave(driveY, onSubSys=0, waves=dragY1(t=20, a=1.7, tau=10, sigma=4))\n",
    "ham.plot()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Simulation and related tools\n",
    "\n",
    "With the Hamiltonian and control pulses defined, we can use the function `simulate()` to calculate the time-ordered evolution operator:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "results = ham.simulate()\n",
    "print(results.result)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "We can also display the detailed information regarding the constructed Hamiltonian using the function `print()` :"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "print(ham)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The pulse sequence for each control term can be extracted by `getPulseSequences()` in `job`:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "print(ham.job.generatePulseSequence(driveY, 0))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Summary\n",
    "\n",
    "The purpose of this tutorial is to introduce how to construct the system Hamiltonian using Quanlse and to simulate and visualize it. After reading this tutorial, users can use this link [tutorial-construct-hamiltonian.ipynb](https://github.com/baidu/Quanlse/blob/main/Tutorial/EN/tutorial-construct-hamiltonian.ipynb) to jump to the GitHub page corresponding to this Jupyter Notebook documentation to get the relevant code, try different parameter values or functions given in this tutorial to gain a deeper understanding."
   ]
  }
 ],
 "metadata": {
  "interpreter": {
   "hash": "2a31ed61199c5c13a03065ecec963a63da8631d96d1b9e695dac4715cb4eadb9"
  },
  "kernelspec": {
   "display_name": "Python 3 (ipykernel)",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.8.12"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 4
}
